rm(list=ls(all=TRUE))

path <- "" #here you need to adjust the path macro to the place on your computer
setwd(path)

library(haven)
library(dplyr)
library(ggplot2)
library(Hmisc)
library(lme4)
library(texreg)
library(effects)
library(estimatr)
library(sjPlot)

LevNum <- function(x){
  c(round(quantile(x, na.rm=T)[2], digits=5),
    round(mean(x, na.rm=T), digits = 5),
    round(quantile(x, na.rm=T)[4], digits=5))
}

dat <- read_dta(paste0(path, "\\party_governments.dta"))

##################################
###### Figure 6, appendix C ######
##################################

dat$iso2code <- as_factor(dat$iso2code)
dat$country <- as_factor(dat$country)

# adjusting the sample, dropping parties which have not sponsored 
# protests and cabinets under which no protest has happened.

table(dat$country)
length(unique(dat$country))

dat <- dat %>%
       group_by(party_id) %>%
       mutate(mean_prts_evts = sum(parties_evts))

dat <- subset(dat, mean_prts_evts!=0)

dat <- dat %>%
  group_by(cabinet_id) %>%
  mutate(mean_prts_evts = sum(parties_evts))

dat <- subset(dat, mean_prts_evts!=0)


# adding some new variables

dat$mainstream[dat$family2==2 | dat$family2==4 | dat$family2==6 | dat$family2==8] <- 1
dat$mainstream[is.na(dat$mainstream) & !is.na(dat$family2)] <- 0
dat$mainstream <- as.factor(dat$mainstream)
table(dat$mainstream, useNA ="always")

dat$past_pm[dat$pm_counter>=1] <- 1
dat$past_pm[dat$pm_counter<1 | is.na(dat$pm_counter)] <- 0
dat$past_pm <- as.factor(dat$past_pm)
table(dat$past_pm, useNA ="always")

dat$opp_party <- as.factor(dat$opp_party)

# low and high contexts

cab_cntr <-  dat %>% 
  group_by(iso2code, cabinet_name) %>%
  summarise_at(vars(all_evts, all_parties_evts), funs(mean(., na.rm=TRUE)))

cab_cntr <- within(cab_cntr, prt_share <- all_parties_evts*100/all_evts)
cab_cntr$cabinet_name[cab_cntr$prt_share>100]
quantile(cab_cntr$prt_share)[3]
mean(cab_cntr$prt_share)


nodiff.cntry <- c("BG", "CY", "CZ", "EE", "DE", "HU", "IS", "LV", "LT", "SK", "CH")
diff.cntry <- c("AT", "BE", "DK", "FI", "FR", "GR", "IE", "IT", "LU", "MT", "NL", "NO", "PL", "PT", "RO", "SI", "ES", "SE", "UK")

cab_cntr$prt_share_cat[as.character(cab_cntr$iso2code) %in% diff.cntry] <- 0
cab_cntr$prt_share_cat[as.character(cab_cntr$iso2code) %in% nodiff.cntry] <- 1

cab_cntr$prt_share_cat <- relevel(as.factor(cab_cntr$prt_share_cat), "0")

dat_upper_lev_cat <- subset(cab_cntr, select=c("iso2code",
                                               "cabinet_name",
                                               "prt_share",
                                               "prt_share_cat"))
dat <- merge(dat, dat_upper_lev_cat)
rm(dat_upper_lev_cat, cab_cntr)

table(dat$prt_share_cat)

# Party level regression

class(dat$prt_share_cat)

dat_diff <- dat[dat$prt_share_cat==0,]
dat_nodiff <- dat[dat$prt_share_cat==1,]

table(dat$all_evts)

m_diff_1 <- lm_robust(parties_evts ~ 
                 state_market + 
                 liberty_authority +
                 mainstream + 
                 past_pm + 
                 opp_party + 
                 vote_bef + 
                 vote_ch_perc_bef +
                 mass_org_scale  +
                 length_gov + 
                 country, data=dat_diff,
                 clusters=country, se_type="stata", alpha = 0.16)

m.eff1 <- tidy(m_diff_1)
m.eff1 <- as.data.frame(m.eff1)
m.eff1 <- subset(m.eff1, select=c("term","estimate",
                                "conf.low", "conf.high"))
m.eff1 <- m.eff1[!grepl("country", m.eff1$term),]
m.eff1$facet <- "differentiated"

m_ndiff_1 <- lm_robust(parties_evts ~ 
                        state_market + 
                        liberty_authority +
                        mainstream + 
                        past_pm + 
                        opp_party + 
                        vote_bef + 
                        vote_ch_perc_bef +
                        mass_org_scale  +
                        length_gov + 
                        country, data=dat_nodiff,
                      clusters=country, se_type="stata", alpha = 0.16)

m.eff2 <- tidy(m_ndiff_1)
m.eff2 <- as.data.frame(m.eff2)
m.eff2 <- subset(m.eff2, select=c("term","estimate",
                                "conf.low", "conf.high"))
m.eff2 <- m.eff2[!grepl("country", m.eff2$term),]
m.eff2$facet <- "not differentiated"

eff.m <- rbind(m.eff1, m.eff2)
eff.m <- eff.m[!grepl("Intercept", eff.m$term),]
eff.m <- eff.m[!grepl("length_gov", eff.m$term),]


eff.m$name[eff.m$term=="state_market"] <- "State-market"
eff.m$name[eff.m$term=="liberty_authority"] <- "Liberty-authority"
eff.m$name[eff.m$term=="mainstream1"] <- "Mainstream"
eff.m$name[eff.m$term=="past_pm1"] <- "Past PM party"
eff.m$name[eff.m$term=="mass_org_scale"] <- "Mass party org."
eff.m$name[eff.m$term=="opp_party1"] <- "Opp. party"
eff.m$name[eff.m$term=="vote_bef"] <- "Vote share"
eff.m$name[eff.m$term=="vote_ch_perc_bef"] <- "Ch. in vote share"

eff.m <- eff.m[!is.na(eff.m$name),]
eff.m$name <- factor(eff.m$name,
                     levels=c("Mainstream",
                              "Past PM party",
                              "Opp. party",
                              "Mass party org.",
                              "Vote share",
                              "Ch. in vote share",
                              "State-market",
                              "Liberty-authority"))

ggplot(data=eff.m, aes(y=estimate, x=name,
                       ymin=conf.low, ymax=conf.high,
                       group=facet,
                       color=facet)) +
  geom_pointrange(aes(shape=facet),
                  position = position_dodge2(width = 0.5, padding = 0.5),
                  show.legend=FALSE) +
  geom_hline(yintercept = 0, linetype="dotted") +
  geom_line(aes(y=conf.low-100)) +
  geom_point(aes(y=conf.low-100, shape=facet), size=2.5) +
  scale_x_discrete(limits = rev(levels(eff.m$name))) +
  coord_flip(ylim=c(min(eff.m$conf.low), max(eff.m$conf.high))) +
  xlab("") + ylab("") +
  theme_bw() +
  theme(legend.title = element_blank())

ggsave(paste0(path, "\\f5_prt_coef_avg.eps"),
       plot=last_plot(), dpi=300, width=7, height=5)
ggsave(paste0(path, "\\f5_prt_coef_avg.png"),
       plot=last_plot(), dpi=300, width=7, height=5)

